pacman::p_load(tidyverse, data.table, sf, raster)
rm(list=ls())
################################################################################ 

###################### Standardize geographic units from shapefiles

shpfolder <- c("../../data/shapefiles/clean/")
country_list <- c("argentina","brazil")

for (country in country_list){
  
  shpfile = sf::read_sf(dsn=path.expand(shpfolder), layer=paste0(country,'_county'))
  df = fread(paste0("../../data/landuse/clean/geographicunits_",country,".csv.gz"))
  
  if(country=='argentina'){
    
    df = rename(df, c("cty_id"="county_id"))
    df = df[,c('cty_id','region','region_id')]
    shpfile_geo = merge(shpfile, df, by='cty_id', all.x=T)
    temp_cw = shpfile_geo[, c("st_id", "region", "region_id")] 
    temp_cw <-`st_geometry<-`(temp_cw, NULL)
    temp_cw = distinct(temp_cw)
    temp_cw = na.omit(temp_cw)
    shpfile_geo = merge(shpfile_geo[,c('cty_id','ctry','st','st_id','cty')], temp_cw, by='st_id', all.x=T)
    
  } else {
    
    df = rename(df, c("cty_id"="county_id", "micro_id"="microregion_id","meso_id"="mesoregion_id"))
    df = df[,c('cty_id','amc_id','micro_id','meso_id','region','region_id')]
    shpfile_geo = merge(shpfile, df, by='cty_id', all.x=T)
    temp_cw = shpfile_geo[, c("st_id", "region", "region_id")] 
    temp_cw <-`st_geometry<-`(temp_cw, NULL)
    temp_cw = distinct(temp_cw)
    temp_cw = na.omit(temp_cw)
    shpfile_geo = merge(shpfile_geo[,c('cty_id','ctry','st','st_id','cty','amc_id','micro_id','meso_id')], temp_cw, by='st_id', all.x=T)
    
  }
  
  st_write(shpfile_geo, dsn = "../../data/shapefiles/clean/", layer=paste0(country,'_all'), driver = "ESRI Shapefile", append=FALSE, quiet=TRUE)

}

#Combine Argentina+Brazil
shpfile_ar = sf::read_sf(dsn=path.expand("../../data/shapefiles/clean/"), layer='argentina_all')
shpfile_br = sf::read_sf(dsn=path.expand("../../data/shapefiles/clean/"), layer='brazil_all')
final_cols = c("ctry","st","st_id","cty","cty_id","region","region_id")
shpfile_ar = shpfile_ar[,final_cols]
shpfile_br = shpfile_br[,final_cols]
shpfile_arbr = rbind(shpfile_ar,shpfile_br)
st_write(shpfile_arbr, dsn = "../../data/shapefiles/clean/", layer='argentinabrazil_all', driver = "ESRI Shapefile", append=FALSE, quiet=TRUE)

#South America
shpfile = sf::read_sf(dsn=path.expand("../../data/shapefiles/clean/"), layer='world_country')
shpfile <- shpfile[shpfile$ctry %in% c('Chile', 'Uruguay','Bolivia', 'Paraguay','Peru','Ecuador','Colombia','Venezuela','Guyana','French Guiana','Suriname'), ]
shpfile$st = shpfile$ctry
shpfile$st_id = shpfile$ctry
shpfile$cty = shpfile$ctry
shpfile$cty_id = shpfile$ctry
shpfile$region = shpfile$ctry
shpfile$region_id = shpfile$ctry
shpfile = shpfile[,final_cols]
shpfile_rosa <- st_crop(shpfile, extent(-85,-32, -65, 11.75)) #bound to mainland South America
shpfile_sa = rbind(shpfile_arbr,shpfile_rosa)
st_write(shpfile_sa, dsn = "../../data/shapefiles/clean/",  layer='southamerica_all', driver = "ESRI Shapefile", append=FALSE, quiet=TRUE)



